Astronomy & Astrophysics manuscript no. dburningarxiv 


©ESO 2012 


October 16, 2012 





Deuterium burning in objects forming 
via the core accretion scenario 

Brown dwarfs or planets? 

P. Molliere 1 and C. Mordasini 1 

Max-Planck-Institut fur Astronomie, Konigstuhl 17, D-691 17 Heidelberg, Germany 
Received 19 June 2012 / Accepted 30 August 2012 

ABSTRACT 

Aims. Our aim is to study deuterium burning in objects forming according to the core accretion scenario in the hot and cold start 
assumption and what minimum deuterium burning mass limit is found for these objects. We also study how the burning process 
influences the structure and luminosity of the objects. Furthermore we want to test and verify our results by comparing them to 
already existing hot start simulations which did not consider, however, the formation process. 

Methods. We present a new method to calculate deuterium burning of objects in a self-consistently coupled model of planet formation 
and evolution. We discuss which theory is used to describe the process of deuterium burning and how it was implemented. 
Results. We find that the objects forming according to a hot start scenario behave approximately in the same way as found in previous 
works of evolutionary calculations, which did not consider the formation. However, for cold start objects one finds that the objects 
expand during deuterium burning instead of being partially stabilized against contraction. In both cases, hot and cold start, the mass 
of the solid core has an influence on the minimum mass limit of deuterium burning. The general position of the mass limit, 13 Mj, 
stays however approximately the same. None of the investigated parameters was able to change this mass limit by more than 0.8 Mj. 
Due to deuterium burning, the luminosity of hot and cold start objects becomes comparable after ~ 200 Myrs. 

Key words, stars: planetary systems - stars: brown dwarfs - planets and satellites: formation - planets and satellites: interiors - 
methods: numerical 



1. Introduction 

It is generally accepted that compact gaseous objects with a 
mass greater than 12-13 Mj, where Mj is the mass of Jupiter, 
will start deuterium burning via the reaction p + d — > y + 3 He 
(Saumon et al. 1996, Chabrier & Baraffe 2000, Burrows et al. 
2001). In masses greater than approximately 63 Mj, lithium 
burning will set in via the two reactions p + 7 Li — > 2a and 
p + 6 Li — > a + 3 He (see e. g. Burrows et al. 2001). As long 
as the object's mass does not exceed a value of approximately 
80 Mj, which is the lower mass limit for hydrogen burning 
(Burrows et al. 2001), it belongs to the class of so-called Brown 
Dwarfs. Objects with masses below 13 Mj, i.e. objects which 
are not able to burn deuterium, are called planets (Boss et al. 
2007), if they are in an orbit around a star or a stellar remnant 
and fulfill the properties which are demanded for planets in 
the solar system as well (independent of the actual formation 
mode). Observationally, for companions in orbit around roughly 
solar like stars, it seems to be difficult to discriminate between 
planets and low-mass Brown Dwarfs: The observed frequencies 
of objects in the mass range around 13 M s behave smoothly 
and exhibit no special pattern which might be expected if 
two different formation processes (one for planets and one for 
Brown Dwarfs) would be at work (Segransan et al 2010). Thus, 
high-mass planets and low-mass Brown Dwarfs might form 
through the same processes. Two different ways of forming 
(giant) planets are being discussed. The first is the so-called disk 
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instability or direct collapse model, which explains the forma- 
tion of planets by a gravitational instability in the protostellar 
disk, whereas the second model is the core accretion model, 
which explains the formation of giant planets with building up 
a solid core by planetesimal accretion which, at some point, 
becomes heavy enough to accrete large amounts of gas from the 
protostellar disk (for a detailed review of both methods see e.g. 
D'Angelo et al. (2010)). In this paper we present a new model 
that combines the effect of deuterium burning with the theory of 
core accretion planet formation (based on the implementation 
by Alibert et al. (2004) and Mordasini et al. 2011). A difference 
to the investigation conducted so far on deuterium burning 
objects is the presence of a solid core in the center of the 
planet. It is interesting to study how the presence of the core 
affects the burning process, as previous work has shown that 
deuterium burning in objects harboring a solid core is possible 
(Baraffe et al. 2008). An advantage of the combined study of 
deuterium burning and planet formation is the the possibility to 
investigate the reaction of deuterium burning to the variation of 
the parameters that constrain the planetary formation process. 
Such parameters can be the dust-to-gas ratio in the protostellar 
disk, the maximum allowed gas accretion rate or the mass 
fraction of helium in the gas etc. As the existence of solid cores 
in giant planets has been subject to a discussion (e.g. Guillot et 
al. (2004) and Wilson & Militzer (2012) suggest that the cores 
might dissolve), our results are for the limiting assumption that 
the core does not dissolve. 

The results of Lubow et al. (1999) suggest that objects 
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forming by core accretion might not get massive enough to burn 
deuterium. The reason for this is the formation of a gap in the 
disk for sufficiently massive objects which prevents them from 
further significant mass accretion. This would limit the maximal 
masses to ~ 10 Mj. However, Kley & Dirksen (2006) find 
that as an object reaches a certain mass threshold, the gap gets 
wide enough to clear the locations of the eccentricity damping 
resonances, such that the resonances exciting eccentricity 
become important. The system is then prone to an eccentric 
instability, resulting in again high accretion rates of the object. 
This eccentric instability is the motivation for considering 
objects forming via the core accretion scenario with masses 
high enough to fuse deuterium. As the exact temporal behavior 
of the mass accretion rate in these phases is still unknown, 
we used the simplification of assuming a constant maximum 
allowed gas accretion rate M max in the runaway accretion phase. 

1.1. Organization of the paper 

The paper is organized as follows: In Section 2 we explain the 
theory which was used to describe the deuterium burning pro- 
cess. In Section 3 we describe how deuterium burning was im- 
plemented into the already existing code. In Section 4 we com- 
pare our results obtained for the deuterium burning assuming a 
hot start to already published work (in which also a hot start was 
assumed, but the formation process was left out). In Section 5 we 
look at the deuterium burning process in objects forming via the 
cold start assumption. In Section 5.2 we investigate the implica- 
tions of changes in quantities such as the dust-to-gas ratio in the 
protoplanetary disk, the initial deuterium abundance, the max- 
imum allowed mass accretion rate or the helium mass fraction 
of the gas for the minimum mass limit for deuterium burning. 
In Section 6 we show how deuterium burning can shorten the 
timescale after which the luminosity of hot and cold start objects 
becomes comparable. Our results are summarized in Section 7, 
followed by the conclusion in Section 8. Finally, analytical ap- 
proximations concerning the cooling phase of the objects can be 
found in Appendix A and B. 



2. Theory of deuterium burning 

In order to describe the deuterium fusion process one has to cal- 
culate its contribution to the luminosity. This means that we have 
to know the energy generation rate e, i.e. energy per unit mass 
and unit time. We use the rate as it is defined in Kippenhahn & 
Weigert (1990), assuming fully ionized, non-degenerate nuclei. 
The value assumed for the energy released in every fusion pro- 
cess was taken from Fowler et al. (1967), whereas the value for 
the cross-section factor was taken from Angulo et al. (1999). 



2.1. Screening 

In order to calculate the nuclear reaction rates of the deuterium 
fusion process correctly, screening theory had to be applied, i.e. 
the effect that enhances the nuclear fusion rate if the repelling 
positive charges of the nuclei are shielded from each other by 
surrounding electrons. We implemented screening following the 
papers of Dewitt et al. (1973) and Graboske et al. (1973). 

2.2. Electron degeneracy 

As we treat objects with high central densities, electron degener- 
acy will become important. The electron degeneracy will affect 



the screening behavior. The factor e accounts for this effect and 
measures the strength of degeneracy. It is defined as (Dewitt et 
al. (1973)): 

1 F_i/#) 

where F v (if/) are the Fermi-Dirac integrals 



now 



dx 



(2) 



and iff is the so-called degeneracy parameter. For if/ — > oo the 
electrons are fully degenerate, whereas for if/ — > -oo they are 
fully non-degenerate. As if/ goes from -oo to oo ® e goes from 1 to 
0. To obtain if/ we tabulated Fx/iiif/) in the i^-range of interest and 
used that it holds that (see e.g. Kippenhahn & Weigert (1990)) 

= - • (3) 
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3. Model 



A detailed description of the simulation method used in this pa- 
per can be found in Mordasini et al. (201 1). Here we only give a 
brief outline of its functioning, concentrating on changes to the 
already existing method. The following four equations needed to 
calculate the object structure hold for spherical symmetric and 
hydrostatic conditions. They are the standard equations for cal- 
culating stellar structures. The assumption of hydrostatic equi- 
librium can be justified as the planet is nearly hydrostatic even 
in the most critical phase of the radial collapse of the envelope. 
This happens at the transition from the attached to the detached 
phase, where one finds the kinetic energy of the gas to be many 
orders of magnitude smaller than the total energy of the object. 
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However, instead of Eq. 5, we use the following equation 

m 

~dr 

and calculate the luminostiy as described in Sect. 3.1. The 
reason why it is justified to choose dl/dr — is outlined in Sect. 
3.1.1. 

Eq. 5 describes the radial evolution of the temperature, 
where the temperature gradient V is defined as 

_ dlnT 



d lnP- (?) 
Following the Schwarzschild criterion and assuming that the ob- 
ject is either only convective or radiative in a shell, we set 

V = min(V ad ,V rad ) (8) 

where V a d, the adiabatic temperature gradient, is given directly 
by the EOS. The radiative temperature gradient V ra d is defined 
as 

V» a = ^^ (9, 



l6nacG mT 4 

In Eq. 9, a denotes the radiation energy density constant, c the 
speed of light, k the opacity and I = l(r) the luminosity. Finally, 
the right part of Eq. 5 yields the luminosity change per unit 
radius, where S denotes the entropy per unit mass, e denotes 
the thermonuclear energy generation rate per unit mass and unit 
time. 
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3.1. Calculation of the luminosity 

Starting from the detached phase 1 , the model, as described in 
Mordasini et al. (201 1), uses a shooting method to find the struc- 
ture of the object, iterating on the outer radius. As one needs 
an surface boundary condition for the luminosity, the following 
derivation is chosen: The total energy of the object is 



F + F- 



r M Gm r A 

dm + 

Jo r Jo 



udm (10) 



where u is the internal energy per unit mass. Motivated by the 
virial theorem and the fact that our solutions are found assuming 
a hydrostatic structure, one can express the total energy of the 
object as 

GM 2 

^ = ~^ (ID 



The factor £ depends on the radial mass distribution inside the 
object and its internal energy content 2 . In order to obtain the lu- 
minosity one can differentiate Eq. 1 1 with respect to time which 
yields 

GM . GM 2 . GM 2 . 
L = £_ M _£ /?+ £ (12) 
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This is only correct, of course, if we do not have the effect of 
deuterium burning, which will be included in the luminosity con- 
sideration in the next subsection. As we never know £ in advance 
of the structure calculation at a certain timestep (and thus also £ 
is unknown), we use the £ of the previous timestep and approxi- 
mate 

/ GM . GM 2 . \ 



R 



2R 2 



The correction factor C accounts for the fact that we excluded 
the £-term from our calculation (which we assume to be small, 
this can be obtained by choosing the timestep sufficiently short) 
and that we used £ from the previous timestep. An estimate of 
the correction factor for the next timestep to t + dt is then found 
by 

_ gtotCj) - E tot (t - dt) 

-L(t)dt ( ' 

where the total energies are obtained via the integrals defined in 
Eq. 10. When considering a hot start, the internal luminosity Lj nt 
(i.e. the luminosity supposed to originate from within the planet) 
used to calculate the internal structure of the planet is simply the 
luminosity as defined in Eq. 13 



In the cold start case we set 
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(15) 



(16) 



assuming that the gravitational potential energy released by the 
gas freely falling onto the planet is radiated away very efficiently 
in a shock at the objects surface. This causes the gas to lose its 
initially quite high entropy and it is incorporated into the object 
at a low specific entropy. One must add, however, that it is still 

1 We start the consideration from the detached phase here as only 
after the collapse and the subsequent runaway accretion of mass one 
reaches temperatures and densities high enough to fuse deuterium. 

2 For an ideal gas with an EOS in the form P = (y - l)pu the de- 
pendence of £ on the internal energy content would be expressed as an 
dependence on y, the adiabatic index of the gas. 



unclear whether objects forming via the core accretion scenario 
really form "cold", as this is depending completely on the struc- 
ture of the shock through which the gas is accreted in the run- 
away accretion phase (Marley et al. 2007, Stahler et al. 1980). 
It might be possible that the gas accretion luminosity should, 
at least partially, be treated as an internal energy source (com- 
monly denoted as a "warm start"), i.e. the gas is accreted with 
a higher specific entropy. Motivated by the fact that the shock 
structure and thus the initial entropy of the freshly accreted gas 
is not really understood, the post-formation evolution of plan- 
ets or objects with initial conditions which lie in the warm start 
regime were studied recently (Spiegel & Burrows (2012)). 

3.1.1. Further simplifications 

As one finds the objects to be fully convective in the detached 
and post-formation evolution phase (except for a thin radiative 
layer at the outer boundary) we do not calculate the radial struc- 
ture of the luminosity, as it is not necessary. The radial structure 
of the luminosity only enters the objects structure calculation via 
Eq. 5 and only if the radial shells under consideration are radia- 
tive (see Eq. 8). As the objects are basically fully convective, 
we set d l/d r - (for simplicity), which is not affecting the 
simulation. Furthermore it has the advantage to also cover the 
objects evolution in the attached phase, were the main luminos- 
ity is due to the solid accretion of the core, thus the luminosity 
is approximately constant in the gaseous envelope. As shown in 
Mordasini et al. (2011) the calculation of the structure as ex- 
plained here leads to evolutionary sequences which are in excel- 
lent agreement with the conventional entropy method (see e.g. 
Burrows et al. 1997, Burrows et al. 2001 and Chabrier & Baraffe 
2000). However, a possible caveat might be that we assume the 
objects to have very effective, large scale convection (i.e. the ra- 
diative losses of rising convective blobs are small compared to 
the energy they transport). This is equal to treating the objects 
as adiabatic in the convective regions. This, however, must not 
necessarily be the case (Leconte & Chabrier 2012), given that 
less effective convection is superadiabatic. 

3.2. Adopting the calculation of the correction factor C to 
deuterium burning 

As outlined in the previous section, the boundary luminosity L mt 
needed for the structure calculation of the objects can be approx- 
imated by using the virial theorem. In this method, however, the 
luminosity contribution of deuterium burning cannot be included 
directly. To overcome this problem we utilized the simplification 
of using the deuterium burning luminosity L D calculated at the 
previous timestep, i.e. at t - dt, and set (for the cold start): 



L int (f) = L(t) ■ 



GM(t)M s . d M 
R(t) 



+ L D (t - dt) 



(17) 



The error of this approximation, similarly to the error we make 
for the approximation of the correction factor C, gets smaller 
the smaller we choose the timesteps. In order to get the correct 
value of the correction factor C one has to take into account that 
the total energy definition in Eq. 10 does not include any term 
considering the deuterium fusion process. Thus, we set 



C 



E tot (t) - E tot (t - dt) 
- (L(t) - L D (t - dt)) dt 



(18) 



in order to correct for this. The fact that Lp(f - dt) appears 
here instead of Lu(t) comes from the just mentioned problem 
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that we cannot "predict" the deuterium burning luminosity with 
the virial theorem ansatz and thus use Ld from the previous 
timestep. C is thus still the ratio of the actual difference of po- 
tential and internal energy and the estimated energy radiated 
away due to contraction, cooling and accretion processes in the 
timestep of length dt. 

3.3. Calculation of the deuterium abundance 

In order to calculate the decrease of the deuterium abundance 
we assumed that the deuterium decreases homogeneously in the 
gaseous envelope. The justification for this is that convective 
mixing is found to be very effective. Using the mixing length 
theory and assuming a mixing length equal to the pressure scale 
height, we found that the convective mixing timescale is always 
at least 2 or 3 orders of magnitudes smaller than the timescale of 
deuterium burning. 

4. Results of the hot start accreting model and 
comparison to hot start simulations 

In order to test the results we obtain for deuterium burning with 
our simulation, especially when varying parameters such as the 
object's helium abundance, deuterium abundance etc., we will 
compare with the results found by Spiegel et al. (2011) (here- 
after called SBM1 1). As SBM1 1 considered hot start initial con- 
ditions for their evolutionary models we will also use the hot 
start condition for the outer boundary condition of the luminos- 
ity as described in subsection 3.1. However, it remains an im- 
portant difference that in our model the formation of the objects 
is taken into account whereas in the model of SBM11 only the 
evolution of the objects is considered. The other outer boundary 
conditions for the evolutionary phase after the formation are 



P = 



[ eqm 
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— 
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(19) 



(20) 



This are the gray photospheric boundary conditions in the so- 
called Eddington approximation. A denotes the albedo of the 
object. The simulation of SBM11 uses a more sophisticated at- 
mospheric model, assuming proper non-gray atmospheres. One 
finds, however, that the results one obtains with our simple treat- 
ment agree well with the results obtained from the non-gray at- 
mospheric model (Bodenheimer et al. (2000)). In Tab. 1 the def- 
inition of our fiducial model is given. Thus, all other runs of the 
simulation will only be specified by the parameters in which they 
deviate from this fiducial model. In Tab. 1, a is the distance to 
the star at which the planet is supposed to form (migration is 
disabled), £ Sj o denotes the initial local surface density of solids, 
r ne b and F ne b stand for the temperature and the pressure in the 
disk at the location of the planet (the disk evolution is neglected) 
and / opa stands for the grain opacity reduction factor, i. e. the 
factor by which we assume the grain opacity of the interstellar 
medium to be reduced (for further details, see Mordasini et al. 
(2012)). The deuterium number fraction [D/H] is defined as the 
ratio of the deuterium and hydrogen atoms. In Tab. 2 the def- 
inition of the runs carried out in order to compare to SBM11 
is given. The definition of the metallicity used in Tab. 2 is the 
following (we assume scaled solar compositions) 



[Fe/H] = log 



/ Dust - to - gas ratio 



in 



0.04 



Table 1. Fiducial model for the in-situ formation and evolution 
calculation 



Quantity / Property 



Value / Status 



a(AU) 


5.2 


S s ,o (g/cm 2 ) 


10 


Mgas,™* (MJyr) 


0.01 


Tncb (K) 


150 


-Pneb (dyn/cm 2 ) 


0.275 


Dust-to-gas ratio 


0.04 


Initial embryo mass (M e ) 


0.1 


Migration 


not included 


Disk evolution 


not included 


Planetesimal ejection 


included 


Simulation duration 


10 10 yrs 


Grain opacity red. factor / opa 


0.003 


Helium mass fraction Y 


0.25 


Deuterium number fraction [D/H] 


2 x 10~ 5 


M» 


1M 



Table 2. Hot start comparison runs 



Run name 


Property 




Y22 


Y = 0.22 




Y24 


Y = 0.24 




Y25 


Y = 0.25 




Y26 


Y = 0.26 




Y28 


Y = 0.28 




Y30 


Y = 0.3 




Dl 


[D/H] = 1 x 10 -5 




D2 


[D/H] = 2 x 10~ 5 




D4 


[D/H] = 4 x 10~ 5 




fpgl 


Dust to gas ratio = 0.04 = [Fe/H] 


= 


fpg2 


Dust to gas ratio = 0.09 = [Fe/H] = 


0.36 


fpg3 


Dust to gas ratio = 0.12 = [Fe/H] = 


0.48 


Al 


Mgas^ax = 10-' MJyr 




A2 


M gas , max = 10~ 2 Me/yr 




A3 


M gas>max = 10~ 3 MJyr 





(21) 



This deviates from the usual definition of the metallicity as the 
protosolar mass fraction is found to be 1.49 % (Lodders 2003). 
In the inner parts of the disk, however, one should have a higher 
dust-to-gas ratio: Due to the different orbital velocities of gas and 
solids the solids will start to drift inwards such that one finds a 
dust-to-gas ratio of approximately 0.04 at 5.2 AU if one con- 
siders a solar like star (Rdzyczka et al. 2004). The choice of 
our fiducial model is mainly motivated by the initial conditions 
Pollack et al. (1996) chose and represents a model which lies in 
the regime of conditions one might expect for an object forming 
via the core accretion scenario. The reason for taking the initial 
deuterium abundance to be 2 x 10~ 5 is motivated by Prodanovic 
et al. (2010), as this is the result they found for the local inter- 
stellar medium. This deuterium abundance of 2 x 10~ 5 has more 
or less become the standard value when investigating deuterium 
burning in Brown Dwarfs. It is also the value used in the fiducial 
model of SBM1 1. The choice we made for the initial conditions 
of the comparison runs where we vary the helium or the deu- 
terium abundance is based on SBM11 as we want to compare 
to their results. Our choice for the metallicity (or dust-to-gas ra- 
tio), however, are made due to a different reason. In SBM1 1 the 
metallicity has an impact on the opacity of the object (as a higher 
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Fig. 1. Temporal evolution of the luminosity (upper left), the radius (upper right), the entropy above the core (lower left) and the 
mass averaged ©^-factor in the envelope (lower right) for hot start objects with masses varying from 10 M } to 30 Mj using our 
fiducial model as defined in Tab. 1. The black dashed line in all plots corresponds to the 13 Mj object. The masses corresponding to 
the individual lines are indicated in the plot. The horizontal green line in the plot of the radius indicates the position of 1 Rj. 



metallicity means more molecules and thus a higher opacity). In 
our model, however, we assume that all the solids accreted by the 
objects are eventually accreted on the core. Thus a higher metal- 
licity (or higher dust-to-gas ratio) will not change the opacity, but 
the core mass of the object. For the opacities we assumed a solar 
composition and used the opacities of Freedman et al. (2008). 
The dust-to-gas ratios we chose correspond to the solar case (i.e. 
[Fe/H] = 0) and to higher metallicities. However, they still lie in 
the feasible regime (Santos et al. 2003, Fischer & Valenti 2005). 
Varying the maximum allowed gas accretion rate is, of course, 
also something we cannot compare to SBM11, as they do not 
consider the formation process. The motivation for studying dif- 
ferent maximum allowed gas accretion rates (M max ) comes from 
the fact that different circumstellar disks will have different prop- 
erties at different stages during their evolution. The disk proper- 
ties such as the viscosity and the surface density determine M max 
(Mordasini et al. 2012). The maximum allowed gas accretion 
rates we chose also lie in the expected regime (Mordasini et al. 
2012). Finally and of importance to qualitatively understand the 
results, is an approximative expression for the nuclear energy 
generation rate (Stahler 1988), 



e oc [D/H] ■ p ■ T ns (22) 



and the following equation, which approximates the central tem- 
perature (see, e. g., Kippenhahn & Weigert 1990): 

M 

Tcent <* — (23) 
K 

This equation is obtained from the equation of hydrostatic equi- 
librium using the EOS of an ideal gas. We found that this expres- 
sion approximates the qualitative behavior of the central tem- 
perature quite well in terms of predicting where it will increase 
during the phase of collapse and runaway accretion. However, 
an important caveat of this approximation is that it assumes an 
ideal gas. This is not true as the objects will become more and 
more degenerate as they contract and their mass increases (see 
Section 5.2 and Appendix B). 

4. 1 . Results of the hot start accreting model: Overall 
behavior 

First we study and compare the overall behavior of the objects, 
forming via hot start core accretion, employing our fiducial 
model. In Fig 1 one can see the luminosity, the radius, the spe- 
cific entropy above the core and the mass averaged degeneracy 
related ® e factor (see Equ. 1). We show the evolution starting 
at 10 6 years, as the foregoing growth of the solid core to the 
isolation mass and the subsequent attached phase of the ob- 
ject are equal in the cold start case of the core accretion scenario. 
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As one can see in the radius plot of Fig. 1, as the planets 
are still growing in mass via runaway mass accretion, they will 
also grow in size (after having collapsed due to the transition 
from the attached to the detached phase) until they have reached 
their final mass. The end of the mass accretion produces 
as sharp bend at approximately 3 x 10 6 years in the radius 
evolution for the 30 Mj object (and earlier for the lower mass 
objects) and inverts the radius growth phase into a contraction 
phase. This is a difference to the cold start model of the core 
accretion scenario, where one observes a gradual contraction 
of the planet once it has collapsed when entering the detached 
phase (see e.g. Mordasini et al. 2011, Marley et al. 2007). This 
difference is due to the inclusion of accretion luminosity into 
the internal luminosity. Thus, the gas is not accreted onto the 
already contracting object as in a cold start, but it is accreted 
injecting all its released gravitational potential energy into the 
object, which prevents it from contraction. This means that the 
gas is accreted at high entropy. The radii one finds directly at 
the end of formation can be regarded as possible initial radii 
for hot start evolution models if one assumes that the objects 
have formed according to the core accretion scenario. However, 
there is one caveat which must be addressed: As the accreted 
gas first deposits it's liberated gravitational energy in the upper 
regions of the object in the runaway accretion phase of a hot 
start, the simplifiying assumptions made for the luminosity in 
Sect. 3.1.1 might break down. During a "hot" accretion radiative 
zones well inside the object could develop (see e.g. Mercer- 
Smith et al. 1984), so that our assumption that d l/d r = 
becomes invalid in these regions. Thus by forcing the object 
to have an adiabatic structure our results could depart quite 
significantly from the correct values in this phase. However, 
as the Kelvin-Helmholtz timescale of an object at the end of 
formation is very short the objects would quickly forget about 
their thermodynamic state directly after the formation, evolving 
back on the correct evolutionary track (provided that the latter 
is not characterized by a much lower entropy). Thus one 
should be cautious when interpreting our results directly after 
formation as initial conditions for a hot start evolutionary model. 

Looking at the specific entropy above the core in Fig. 1 
one can see how the most massive object (30 Mj) reaches a 
specific entropy as high as 14.5 £ B /Baryon due to the accretion 
of high entropy material (which is characteristic for a hot start). 
We note here as well that these values found for the specific 
entropy at the end of formation might be realistic initial values 
for evolution calculations if one assumes that the objects have 
formed via core accretion. However one should keep in mind 
the caveat outlined above. As the objects are fully convective in 
their interior (and thus have an adiabatic temperature gradient) 
the specific entropy is constant throughout the envelope. Again 
one sees that also the decrease of the entropy is partially 
stabilized by deuterium burning. 

In the & e plot of Fig. 1 one sees how the degeneracy re- 
lated factor ® e (see Eq. 1) evolves. For this plot & e was mass 
averaged in the gaseous envelope. As one can see the envelope 
starts fully non-degenerate (0 C = 1) and reaches a first minimum 
at approximately e = 0.7. Instead of starting at 10 6 years as all 
other plots of Fig. 1, the e plot starts at 10 5 years in order to 
show the initial value of f = 1 . After having reached a value 
of 0.7, & e then starts to increase again at the beginning of the 
runaway phase and then drops down to 0.6 at the moment where 
the objects envelope collapses (at « 10 6 years). One can see that 



deuterium burning will partially stabilize & e as well. After this 
stabilization phase, however, & e will start do decrease stronger 
again, approaching small values between and 0. 1 at the end of 
the simulation, corresponding to stronger degeneracy. 



4.1 .1 . Comparison of the stabilization radii to the results of 
Burrows et al. (1997) 

As the SBM11 hot start calculations (to which we compare our 
work in the next section) are based on the model of Burrows et al. 
(1997), it seemed worthwhile to compare to those results, espe- 
cially as the data of the Burrows et al. (1997) is available online. 
In Fig. 2 one can see the radius evolution obtained by Burrows 
et al. (1997) (dashed green lines) and our own results (red solid 
lines). The results of Burrows et al. (1997) have been shifted in 
time in order to make the phases of deuterium burning coincide. 
This is justified as the Burrows et al. (1997) results do not con- 
sider the formation of the object, thus the definition of t = is 
not equal. As one can see in Fig. 2, the radius at which the par- 
tial stabilization of the objects occurs is the same in both models 
and also the subsequent radial evolutions are nearly equal. This 
again underlines the reliability of our model assumptions, i.e. the 
way how we calculate the luminosity and the way how deuterium 
burning was included into the simulation. 

4.1.2. Comparison to the SBM11 results 

Comparison with SBM1 1 (SBM1 1 also assumed a Helium mass 
fraction of 0.25, along with a deuterium number fraction with 
respect to hydrogen of 2 x 10~ 5 as their fiducial model) shows 
that the radii at which the planets become partially supported 
against contraction by deuterium burning roughly coincide with 
our results. A 20 Mj object will be in this phase at a radius 
between 2 and 2.5 Jupiter radii in the SBM11 results, which 
agrees with our findings. For the lower mass objects in the 14 
to 15 Mj range this phase lies around 1.5 Jupiter radii which is 
also true for our results. 

Furthermore, when looking at the evolution of the ratio 
Ld/L <o1 we find that for a 13 Mj object it peaks at approximately 
0.5 which is the same result SBM11 obtain. Also, we and 
SBM11 both find that for M > 14 Mj the peak of L D /L tot is 
above 0.8, reaching more than 0.9 for the heavier objects > 20 
Mj. The results for L D /L tot for objects below 12 Mj are also 
similar to the SBM11 results, the only deviation is observed 
for the 12 M j object, where our peak value for Lo/itot is twice 
as high. One must note, however, that L D /L tot is a very steep 
function of the mass in this range, thus slight differences in the 
models might be seen quite strongly. 

Also, when looking at the luminosities at which the ob- 
jects are partially stabilized due to deuterium burning we find 
a very good agreement between our results and those of SBM1 1 . 

In our model and in the model of SBM11 the metallicity 
enters in two different ways: It sets the opacity in the SBM11 
case and the core mass in our case. As both SBM11 and we 
ourselves find that the metallicity has a significant influence 
on the deuterium burning process (see Section 4.2), it is rather 
a coincidence that we straightaway found the metallicity for 
our fiducial model such that it yields approximately the same 
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Fig. 2. Radial evolution of deuterium burning objects (hot start accretion). The solid red lines show our own results, whereas the 
dashed green lines show the results of Burrows et al. (1997). As the t — definition is not the same, the Burrows et al. (1997) results 
where shifted in time in order to make the phases of deuterium burning coincide. The mass of the objects is indicated in the plots. 



deuterium burning behavior as in the SBM11 case. Thus 
even though they also assumed a solar metallicity in their 
fiducial model it was per se not clear that we would get a good 
agreement. It shows, however, that once one has found this 
metallicity value which seems to correct for the effect of the 
different implications of the metallicities in both models, the 
results seem to be approximately the same in all considered 
quantities. 

4.2. Parameter study 

In their paper, SBM1 1 tested how the deuterium burning process 
is influenced when parameters such as the hydrogen abundance, 
the deuterium abundance and the metallicity are varied. We also 
carry out this parameter study to further test our own simulation, 
before we proceed to the study of the cold start results of 
objects forming by core accretion in Section 5. In the upper 
left and right panel of Fig. 3 one can see the Lo/^tot evolution 
and the deuterium depletion for the runs with different helium 
abundances (Y22-Y30) for a 13 Mj object. 

The bumpy pattern in some Lo/L tol plots is due to the in- 
terpolation of the opacities in the tables of Freedman et al. 
(2008) and thus a numerical artifact. 

Concerning the height of the Lu/L tol peaks we again find 
a very good agreement with the results of SBM11. We do, 
however, also observe that our objects burn 10-50 % more 
deuterium as they sustain a high burning rate over a longer 
period of time. Overall, we thus confirm the fact that a higher 
helium abundance yields a stronger deuterium burning process, 
due to the increased central density of the objects (see Eq. 22) 3 . 
The results of the runs with different deuterium abundances 
(Dl, D2 and D4) are shown in middle row of Fig. 3. Again, we 
find a very good agreement when comparing to the peak value 
of Lu/L tot of the Spiegel et al. (2011) results. Furthermore, 
we again find the already mentioned fact that our objects burn 
more deuterium than in the SBM11 case. A possible reason for 
this could be the use of different opacities: Higher opacities 
enable the objects to keep high central temperatures over a 

3 In the terminology used in this paper, by saying "central" we always 
mean the innermost layer of gas just above the solid core 



longer period of time (the so-called blanket effect). Furthermore 
SBM11 use a proper atmospheric model while we use gray 
photospheric boundary conditions. As the deuterium burning 
rate is depending on the temperature very strongly (see Eq. 22) 
this can sustain the burning rate over a longer period of time 
as well. We conclude that we are in a good agreement with 
the results by SBM11 but out objects bum up to 50 % more 
deuterium. Finally we consider the runs with varied dust-to-gas 
ratio (i.e. varied metallicity) and with varied maximum allowed 
gas accretion rates, which we cannot compare to the SBM11 
results, as they do not consider the formation process or the 
presence of solid cores. It is worthwhile to study the impact of 
these parameters nonetheless. In the lower row of Fig. 3 one 
can see the results obtained for the L D /L tot evolution and for 
the deuterium depletion for runs with different dust to gas ratios 
fpgl, fpg2 and fpg3, corresponding to dust-to-gas ratios of 0.04, 
0.09 and 0.12, respectively. The considered dust-to-gas ratios 
have a significant influence on the deuterium burning behavior. 
The higher it is, the higher the nuclear burning luminosity 
compared to the total luminosity and the more deuterium is 
burned. One can, in principle, imagine at least two reasons for 
this behavior. 

One reason is that a higher dust-to-gas ratio will lead to a 
higher solid surface density in the disk, which will in turn lead 
to a higher solid accretion rate. A higher solid accretion rate 
also means a higher core luminosity, which might influence the 
objects structure as the core accretion luminosity is included into 
the structure calculation. This is motivated by the assumption of 
our model that the solids impacting the object get destroyed and 
the debris then sinks to the core, following the so-called sinking 
approximation (Pollack et al. 1996). Thus one would expect the 
central temperature to be higher, which has a large influence 
on deuterium burning, as the nuclear energy generation rate is 
highly dependent on the temperature (see Eq. 22). 

The other possible reason is due to the fact that we as- 
sume that all solids accreted by the object eventually get 
incorporated in the core, leading to a higher core mass. A higher 
core mass will, even though the core is also growing in size as 
it grows in mass, result in a higher gravitational acceleration on 
top of the core. In order to remain in hydrostatic equilibrium one 
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Fig. 3. Left column: The ratio Lu/L tot . Right column: The depletion of deuterium. The first row shows our results for different 
helium abundances, the second row shows the influence of different deuterium abundances and the last row shows different dust-to- 
gas ratios (as defined in Tab. 2). The mass of all objects is 13 M s (hot start). The run names corresponding to the individual lines 
are indicated in the plot. 



thus needs a higher pressure gradient, i. e. a higher temperature 
gradient, as the temperature gradient in a convective layer 
(which the central gas layer surely is) is 



d T 
37 



1 d lnp Gm(r) 
cp d lnr r 2 



(24) 



which is proportional to the gravitational acceleration. This 
means that a higher temperature gradient might lead to a higher 
central temperature, if it dominates over the fact that the gas 
cannot reach as deep within the object as it could with a lighter 



core (as objects with a higher core mass are found to have a 
smaller total radius and a larger core radius). If it does, this 
would create a higher deuterium burning rate. 

We were able to identify that it is rather the increased core mass 
than the additional core accretion luminosity which affects the 
deuterium burning process in the following way: In Fig. 4, one 
can see the Lo/L tot evolution and the deuterium depletion for 
the runs with a different maximum gas accretion rate (A1-A3). 
As one can see, the object with the lowest maximum mass 
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runs (right plot) (hot start). 



accretion rate reaches the the peak value the latest, as it needs a 
longer time to form. However, all three runs have the same peak 
value of L D /L tot , and the same subsequent evolution of L D /L tot 
after having reached the peak value. One can also see that all 
three runs burn approximately the same amount of deuterium. 
The remaining slight variations in the A1-A3 runs for the peak 
values of Lo/L tot and the deuterium depletion are mainly due 
to a simulation artifact: As we cannot switch off the accretion 
instantaneously (in order to remain numerically stable) it is 
difficult to exactly reach the desired total mass. The deuterium 
burning rate depends on the objects density (see Eq. 22) and 
thus one can see that a slightly different final mass will have an 
impact. 

In conclusion we see that the amount of deuterium burned 
and the strength of deuterium burning is independent of the 
luminosity due to gas accretion, which is included in the 
internal luminosity in the hot start case (see Eq. 15). We have 
furthermore found that the contribution of the solid accretion 
luminosity to the total accretion luminosity is very small in the 
phase of runaway gas accretion. This makes it very unlikely 
that the different solid accretion luminosities in the runs with 
different dust-to-gas ratios ("fpg"-runs) have a strong impact 



on the deuterium burning process. The luminosity due to gas 
accretion is much stronger, but apparently has very little effect 
on deuterium burning. It rather seems to be the presence of the 
core itself which influences deuterium burning via Eq. 24. 

In order to further show this, we plot in Fig. 5 the central 
temperature for the different dust-to-gas ratios and for the 
different maximum gas accretion rates. One clearly sees that for 
the runs with different dust-to-gas ratios, the central temperature 
is always higher for the runs with higher dust-to-gas ratio, 
which persists until the end of the simulation at 10 10 years. 
The runs with different maximum gas accretion rate, however, 
reach approximately the same central temperature at the end of 
formation and their subsequent central temperature is the same. 
The core masses for the three runs "fpgl", "fpg2" and "fpg3" 
are 35.2 M & , 77.5 M e and 101.5 M @ , respectively. For the three 
runs "Al", "A2" and "A3" the core masses are 39.1 M e , 35.2 
M ffl and 31.2 M @ . 



5. Cold start scenario 

In this main section, we consider the formation of objects form- 
ing via the core accretion scenario using the cold start assump- 
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Fig. 6. Formation and subsequent evolution of a 16 Mj object forming via core accretion with a radiatively efficient shock (cold 
start). The four panels show the temporal evolution of the total luminosity (internal + shock) (upper left), the radius (upper right), 
the central temperature (lower left) and the surface temperature (lower right). The dotted blue line shows the results obtained for the 
fiducial model as specified in table 1 . For the solid red line, the deuterium abundance was kept constant, while for the dashed green 
line, deuterium burning was disabled. 



tion, which means that all gravitational energy liberated at the 
accretion shock is radiated away, not contributing to the internal 
luminosity (see Eq. 16). This is equivalent to the accretion of 
low entropy material. We stress again, however, that it is still an 
open question whether the total accretional luminosity is radi- 
ated away at the shock (Marley et al. 2007, Stahler et al. 1980). 
A cold start is the standart scenario for the formation via core ac- 
cretion with a gradual accretion of matter. The typical timescale 
for the runaway gas accretion phase during the formation of a 
13 Mj object is approximately 4 x 10 5 years if one considers a 
M max of 10~ 2 M e /yr (almost the complete mass of the object is 
accreted in the runaway phase). This is quite long compared to 
the typical formation timescale of direct collapse models, which 
is approximately equal to the orbital timescale and of the order 
of 1 2 - 1 3 years in the distances feasible for a direct collapse to 
occur (fldkect > 50 AU for 10 Mj) (Janson et al. 201 1). It is, how- 
ever, unclear if a clump at final mass forms in the direct collapse 
model or whether it would accrete gas as well. For our fiducial 
model, we will use the parameters specified in Tab. 1 . In Section 
5.1, we study the general behavior of objects forming with a cold 
start. In Section 5.2 we investigate the minimum mass limit for 
deuterium burning. In order to do so we vary the parameters of 
the model and carried out the runs given in Tab. 3. 



5.1. Overall behavior 



In order to understand the implications of deuterium burning on 
the evolution of objects forming via core accretion and a radia- 
tively effective shock we first describe the qualitative behavior 
studying the formation and subsequent evolution of a 16 Mj 
object. In order to get an understanding for the changes deu- 
terium burning introduces to the internal structure we also con- 
sider the artificial cases where either deuterium burning is dis- 
abled, or where the deuterium abundance is kept constant. Fig. 6 
shows the evolution of the object's luminosity, it's radius and it's 
central and surface temperature for the three cases (deuterium 
burning, no deuterium burning, constant deuterium abundance). 
Differences in the evolution can only be seen after the object's 
formation is completed, as marked by the end of runaway gas 
accretion which produces the sharp decrease in the luminosity 
at approximately 1.3 x 10 6 years. This is due to the disappear- 
ance of the accretion luminosity. It is clear, however, that the 
phases preceding the gas runaway accretion cannot be important 
for deuterium burning, as the temperatures as well as the densi- 
ties in the gaseous layers are far too low. 
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5.1.1. Luminosity 

Concentrating on the panel with the luminosity in Fig. 6, we see 
that after the formation of the object, the deuterium burning will 
cause the object to behave quite differently. We first note that 
deuterium burning produces an additional, rather flat peak in lu- 
minosity at approximately 10 7 years where the luminosity is an 
order of magnitude higher than in the case with deuterium burn- 
ing disabled. We also see that for objects where the deuterium 
abundance is kept constant, the objects will reach a (hypotheti- 
cal) deuterium burning main sequence, where the object's lumi- 
nosity is generated by deuterium burning only. This is what one 
would naively expect when keeping the deuterium abundance 
constant. 

5.1.2. Radius 

Comparing to the hot start scenario (see Fig. 1), a clearly differ- 
ent behavior can be seen for the evolution of the radius. Instead 
of stabilizing the objects against further contraction as in the hot 
start case, the onset of deuterium burning is so strong that the 
object is inflated from 1.2 to approximately 1.7 R s for the nomi- 
nal case. This suggests that a large part of the deuterium burning 
luminosity is converted into gravitational potential energy in that 
phase, which corresponds to a L D /ij n t ratio (much) larger than 
1 (we will discuss this later). Thus, instead of contracting and 
cooling at high temperatures and gradually burning deuterium 



already at quite large radii as in the hot start case, a cold start 
yields a more abrupt start of deuterium burning. This is due to 
the effect that the mass is accreted at low entropy onto an object 
with high density, due to the collapse from the attached to the 
detached phase. Thus the objects have then much higher burn- 
ing rates as the objects are already very dense when they reach 
the temperature needed for deuterium burning (see Eq. 22). The 
expansion process will be discussed in more detail below. 

5.1.3. Temperatures 

Finally, one can see that deuterium burning also produces a flat 
peak in the central and in the surface temperature when com- 
pared to the case without deuterium burning. The case where the 
deuterium abundance was set constant yields a main sequence- 
like behavior again. 

5.1.4. Expansion 

In order to understand the expansion we compare the radii to 
which the objects expand in the cold start scenario with the radii 
at which the objects get partially stabilized against contraction 
in the hot start case. A comparison like this seems sensible as 
it is well known that nuclear fusion processes in stars, such as 
hydrogen burning, have a thermostatic behavior (also known as 
the pressure temperature thermostat or stellar thermostat) (see, 
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Table 3. Parameter study for cold start simulations 
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e.g. Palla et al. 2002). This can also be seen very clearly in the 
runs with constant deuterium abundance, as the objects, driven 
by deuterium burning, attain a stable state in the hypothetical 
deuterium main sequence. Given the thermostatic nature of deu- 
terium burning, we expect the cold start objects to expand to 
approximately the same radii found for the hot start objects in 
the phase where they get partially stabilized. More precisely, as 
the objects in the hot start case are unable to completely stabi- 
lize themselves against further contraction (as Lu/L tot is always 
smaller than 1) we expect that the cold start objects will not be 
able to expand to the same radii as the "quasi-stabilization"-radii 
of the hot start objects, but to a somewhat smaller radius. The 
reason for L£,/L tot < 1 in the hot start case, as well as for the 
somewhat smaller expansion radius in the cold start case is the 
decreasing deuterium abundance inside the objects. When look- 
ing at Fig. 7 we can see this behavior: The cold start objects 
expand to radii somewhat smaller than the radii at which the hot 
start objects are partially stabilized. When this radius increases 
for the hot start objects (i.e. for higher masses) it also does so 
for the cold start objects. For the 16 M s object in Fig. 7, we have 
also plotted the radius evolution for hot and cold start with deu- 
terium abundances held artificially constant. As one can see, they 
both get stabilized at the same radius (which we call "main se- 
quence radius" in the plot). The hot start object contracts to this 
radius from above, whereas the cold start object is at a smaller ra- 
dius when deuterium fusion ignites and must thus expand to the 
main sequence radius. One can observe that the main sequence 
radius approximately coincides with the "quasi-stabilization" ra- 
dius of the hot start object. Thus, the consideration of the in prin- 
ciple unphysical (artificially held constant deuterium abundance) 
main sequence radius and the thermostatic nature of nuclear fu- 
sion processes help to understand not only the stabilization be- 
havior of the hot start objects, but also the expansion behavior of 
the cold start objects. 



5.1 .5. Evolutionary sequences for objects formed by core 
accretion, assuming a cold start 

Similar to Fig. 1 for the hot start case we also show the temporal 
evolution of the luminosity, the surface temperature, the radius, 
the gravitational acceleration at the surface, the specific entropy 
and the mass averaged degeneracy parameter ® e for the cold 
start case in Fig. 8. Except for ® e , the evolution is shown 
starting at 10 6 years, i.e. shortly before the final mass is reached, 
and therefore shows the final stages of the object's formation 
in the runaway accretion phase. The time axis of ® e starts at 
10 5 years in order to see that the envelope indeed was fully 
non-degenerate in the beginning. 

The luminosity shown in the plot is the total luminosity, 
i.e. the sum of the internal luminosity and the accretion lumi- 
nosity, as this is the luminosity which would be seen by an 
observer. The sharp drop in luminosity (from approximately 0. 1 
L Q for the 23 Mj object) marks the end of runaway gas accretion 
and thus the end of the formation process. One can see that for 
the lower masses up to 13 Mj there is no additional bump in the 
luminosity. Instead it decreases monotonously. For the objects 
of higher masses, however, the luminosity does not decrease 
monotonously but has an additional bump due to deuterium 
burning. For almost all masses shown in the plot the maximum 
luminosity due to deuterium burning is reached far after the 
end of formation. For the 23 Mj object this is not the case: The 
luminosity does not increase strongly anymore after the end 
of formation and indeed one finds that the 23 M s object burns 
almost all of it's deuterium during the formation phase (this is 
discussed in more detail in Section 5.1.6). One can also observe 
that the higher the mass of the object the higher is the bump in 
luminosity. 

A similar behavior is found for the surface temperature: 
The objects of masses up to 13 Mj do not exhibit an additional 
bump due to deuterium burning. All objects with masses 
above, however, do. The peak value of the increased surface 
temperature is increasing with mass. Directly after formation 
the 23 Mj object is over 1000 K hotter than the 13 Mj ob- 
ject. As well as in the luminosity case the peak value of the 
deuterium burning induced temperature bump is reached later 
for the lower mass objects. The time difference between the 
23 Mj object and the 14 Mj object is approximately 5 x 10 7 years. 

Everything which was said about the luminosity and the 
surface temperature also applies for the radius, the entropy and 
the degeneracy factor ® e . Deuterium burning introduces an 
additional bump in the temporal evolution of these quantities 
for masses above 13 M } . The peak value of these bumps will 
be reached later for the less massive objects and the peak value 
will be smaller. 

For the surface gravity it also holds that no deuterium burning 
caused effect can be seen for masses up to 13 Mj. For the more 
massive objects deuterium burning produces a dip, rather than 
a bump, as the gravitational acceleration is proportional to R~ 2 , 
where R is the radius of the object. 

Finally it is important to note that the fact that no addi- 
tional bump or dip is seen for masses up to 13 Mj does not 
necessarily mean that these objects do not burn deuterium at 
all. As we will discuss in the following sections, objects with 
masses around 13 Mj already burn significant amounts of their 
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Fig. 8. Temporal evolution (for cold start objects) of the luminosity as seen by an observer, i.e. including the accretion luminosity 
(first row, left), the surface temperature (first row, right), the radius (second row, left), the gravitational acceleration at the surface 
(second row, right), the specific entropy above the core (third row, left) and the degeneracy related factor & e (third row, right) for 
objects with masses varying from 10 Mj to 23 Mj using our fiducial model as defined in Tab. 1. The black dashed line in all plots 
corresponds to the 13 Mj object. The horizontal green line in the plot of the radius indicates the position of 1 Rj. The masses 
corresponding to the lines are indicated in the plot. 



deuterium. The effect of deuterium burning does not introduce 
an additional bump or dip but it slows down the monotonous 
decrease (or increase for log (g)) of the discussed quantities. 



5.1 .6. Implications of the gas runaway accretion phase 

We now concentrate on the behavior of deuterium burning and 
how it is influenced by the mode of formation of the objects, 
which is a core accretion formation with a radiatively efficient 
shock (cold start). As we have already seen, a major difference 
comes from the fact that cold start objects start deuterium 
burning at a much smaller radius than in the hot start case. 
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masses as indicated in the plot. The blue solid lines show the 
total luminosities of the objects (i.e. the luminosity an observer 
would see), whereas the dashed red lines denote the luminosity 
due to deuterium burning only. The luminosities for the different 
masses have been separated in the y-axis by an arbitrary offset 
for better visibility. Up to approximately 1.3 x 10 6 years, where 
the gas runaway accretion ends for the 13 Mj object, all objects 
have the same luminosity. 



Then they expand approximately to the radius they attain in a 
hypothetical deuterium main sequence, due to the thermostatic 
nature of deuterium burning. In Fig. 9 one sees the evolution of 
the total luminosity L tot and of the luminosity due to deuterium 
burning (Ld) for objects of different masses. When studying 
the plot one can notice that the end of gas runaway accretion 
(marked by the steep decrease in the total luminosity of the 
objects) causes a bend in Ld for all masses, which is clearer 
for the lower masses. For the lower masses up to 17 M s , the 
growth of Ld, which is strongly increasing in the gas runaway 
phase, is brought to a halt and remains approximately constant. 
For the heavier masses, L D will decrease after the end of the 
gas runaway accretion. Furthermore we observe that for the 
higher masses (« 23 Mj), Ld reaches a peak still during the 
phase of gas runaway accretion, then decreases slightly before 
showing the bend at the end of gas runaway accretion, after 
which it decreases even faster. This behavior can be explained 
with the equations 22 and 23. After the object has collapsed at 
the transition from the attached to the detached phase, it is not 
yet massive enough to bum deuterium. Furthermore, due to the 
cold start assumption the (major) part of the planets mass is then 
added onto the planet at low entropy, meaning that it grows in 
mass while contracting even further. As the central temperature 
is proportional to M/R (Eq. 23) and the deuterium burning rate 
is depends on the temperature very strongly and linearly on 
the density (Eq. 22), this process of adding mass on the object 
will cause the observed increase of L D in the phase of runaway 
mass accretion. The sharp bend in Ld at the end of gas runaway 
accretion nicely shows the dependence of Ld on the process of 
runaway accretion (as the object has stopped mass accretion L D 
ceases it's strong increase). For the objects of heavier masses, 
which exhibit a peak of deuterium burning inside the phase of 
gas runaway accretion, another effect becomes important. As 
one can see in Eq. 22, the fusion rate of deuterium burning 
is proportional to the deuterium abundance. Indeed we find 



that those objects bum most of their deuterium already during 
their formation phase (as higher masses mean a higher central 
temperature and thus a higher deuterium burning rate). Thus, 
at some point, Ld will not increase further even though the 
central temperature is still rising, as the effect of the decreasing 
deuterium abundance becomes stronger. The reason for the 
following modest decrease of Ld in the runaway phase is that 
the still increasing central temperature (due to Eq. 23) will 
partially counterbalance the impact of the decreasing deuterium 
abundance. 

Finally, as already said above, the objects expand due to 
deuterium burning. Thus Ld must overshoot the intrinsic 
luminosity (to which it contributes) in the expansion phase. 

5.2. Sensitivity of the minimum mass limit of deuterium 
burning 

Similar to the paper of SBM11, we will investigate the im- 
plication of certain quantities on the process of deuterium 
burning, namely in terms of how they affect the minimum 
mass limit for deuterium fusion. In order to do so, we assume 
a core accretion formation, use the cold start assumption and 
consider the runs as specified in Tab. 3. As SBM11 pointed 
out, no exact border can be drawn between deuterium burning 
and non-deuterium burning objects, as the amount of burned 
deuterium varies smoothly as a function of mass instead of 
exhibiting a sharp transition. In this paper, if an object decreases 
it's initial deuterium abundance by more than 50 %, we will call 
it for simplicity a deuterium burning object and objects below 
this border will be called non-deuterium burning objects. The 
mass at which this transition happens will be called M 5( ). One 
should however bear in mind that this choice of the position of 
the border is partially arbitrary, even though 50 % seems to be a 
sensible choice, as it is in the middle of the two cases of (almost) 
no deuterium burning and complete deuterium combustion. 
This illustrates the problem which arises if one persists on the 
rigorous distinction of planets and Brown Dwarfs in the mass 
regime considered henceforth. 

In Fig. 10 one can see the results we obtain for the runs 
in Tab. 3. The fraction of the initial deuterium which has been 
burned until the end of the simulation (10 10 years) is plotted 
against the total mass of the objects. As one can see, the 
transition from objects which burn no deuterium to objects 
which burn all their deuterium has a width of approximately 5 
Mj and takes place roughly between 10 and 15 M } . 

As one can see in Fig. 10, the largest impact on changing 
M 5 o has the variation of the helium abundance. It must be 
noted, however, that a helium abundance as low as 0.22 is 
rather improbable, as the helium fraction due to the primordial 
nucleosynthesis was approximately 0.25 (Spergel et al. 2007) 
and has, since then, increased in time due to hydrogen fusion 
in stars. The top and bottom right panel of Fig 10 show that 
a variation of the metallicity, as well as the variation of the 
deuterium abundance, will shift the place of the transition region 
between deuterium burning and non deuterium burning objects. 

When considering the maximum mass accretion rate M max 
(bottom left panel) we find that, opposite to the results obtained 
for the hot start case, M max now also has an impact on the pro- 
cess of deuterium burning. The reason for this is the following: 
The objects collapse at the transition from the attached to the 
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Fig. 10. Fraction of the initial deuterium burned as a function of mass until the end of the simulation (10 10 years) for cold start 
objects. Top left: Results for different helium fractions. Top right: Results with different dust-to-gas ratios (or metallicities). Bottom 
left: Results for different maximum gas accretion rates. Bottom right: Results for different deuterium abundances. The uppermost 
run names in the plots correspond to the leftmost curves, the next lower run names correspond to the second leftmost curves etc. 
For the Y30 run in the upper left panel the positions of the individual mass runs are indicated with filled brown circles. For the sake 
of clarity we omitted plotting the individual mass run positions for all other runs. Instead a spline function was plotted through the 
individual mass run positions. 



detached phase (as can be seen at approximately 10 6 years in the 
upper right plot in Fig. 6 showing the radius evolution). After 
this collapse the objects will contract further while they accrete 
most of their mass during the runaway accretion process. This 
is in contrast to the hot start case, where the objects expand 
again during the runaway accretion. It is clear that the higher the 
maximum allowed mass accretion rate, the shorter is the time 
the object will have to contract during the runaway accretion 
process. The consequence of this is the following: At higher 
M mdx , the gas accretes onto a larger object, which means that 
the free fall velocity of the gas at the moment it reaches the 
surface of the object will be smaller, as the free fall velocity for 
the infalling matter is 



2GM 



R 



(25) 



M denotes the mass of the object, R is the radius of the object and 
G the gravitational constant. The lower free fall velocity would 
result in a smaller accretion luminosity if M max was constant (c.f. 
Bodenheimer et al. 2000) 



^acc — M max k^p 



vl = 



GMM n 
R~ 



(26) 



The lower accretion luminosity per accreted mass unit means 
that less liberated gravitational potential energy is radiated away 
at the shock, or in other words, that matter of higher entropy is 
incorporated into the object. This can be understood by the fact 
that although the object forms according to a cold start the ma- 
terial is incorporated into the object at larger radii. Thus it must 
have a larger specific entropy, just as the initial radius of hot 
start simulations is increasing with the chosen initial entropy for 
fully adiabatic objects. We find that the accretion of matter with 
a higher entropy means that the object will have a higher cen- 
tral temperature, leading to higher deuterium burning rates. The 
effect of a higher central temperature due to a higher specific 
entropy might surprise in view of the classical result obtained 
for the central temperature of a self-gravitating sphere of ideal 
gas (see Eq. 23): A higher specific entropy in fully convective 
(i.e. adiabatic) objects generally implies a larger radius (see also 
Spiegel & Burrows 2012). This, in turn, should cause a lower 
central temperature (see Eq. 23). In the cases considered here, 
however, partial degeneracy already plays an important role. The 
fact that the degeneracy becomes stronger for smaller radii as 
well as the fact that fully degenrate objects cool as they shrink 
(e.g. white dwarfs) underlines that a larger radii due to a higher 
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Table 4. Fitting parameters for M50 
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Fig. 11. Central temperature of a 7 M } planet as a function of 
it's specific entropy showing the complete post-formation evolu- 
tionary phase. 



specific entropy must not generally imply a smaller central tem- 
perature. Also the presence of the core can cause the envelope to 
cool. Both arguments are outlined analytically in the appendix 
A and B. In Fig. 1 1 one can see see how the central tempera- 
ture decreases as the entropy decreases in a 7 Mj planet in the 
post-formation evolutionary phase. 

5.2.1 . Detailed study of the behavior of M 50 

Fig. 12 shows how M 5 q varies due to changes in the helium 
abundance Y, the maximum allowed runaway gas accretion rate 
M mdx , the deuterium abundance [D/H] and the metallicity [Fe/H] 
(corresponding to different dust-to-gas ratios). The results ob- 
tained for cold start objects are plotted in blue, while the results 
obtained for hot start objects are plotted in red (we varied the 
fiducial model for the hot start objects in the same way as we 
varied it for the cold start objects given in Tab. 3, but suffixed the 
run names with an "h" instead of a "c"). For Y, [D/H] and [Fe/H] 
the results were fitted to a line using the least square method: 



M 50 /M] = m Y ■ Y + b Y 

M 50 /Mj = m [D/H] ■ [D/H] + b [D/H] 

M 50 /Mj = m [Fe / H ] ■ [Fe/H] + % e/H] 



(27) 



For M raax (note the logarithmic x-axis in Fig. 12 for M max ) the 



results were fitted to the function 



M 50 /Mj = logi \cm 



1 M ffi /yr 



(28) 



This means that in a first approximation M50 depends linearly 
on the helium abundance Y, the deuterium number fraction 
[D/H] and the metallicity [Fe/H] . For M max and the dust-to-gas 
ratio (see Eq. 21), however, there is a logarithmic dependence. 

The upper right panel of Fig. 10 shows the metallicity de- 
pendence of M50. However, instead of labeling the x-axis with 
"[Fe/H]" it was labeled with "[Fe/H]] oc gas ". The intention of this 
is the following: In our model the metallicity sets the dust-to-gas 
ratio via Eq. 21, which is only valid at approximately 5.2 AU 
(see discussion of Equ. 21). The dust-to-gas ratio (fpg), in turn, 
sets the local initial solid surface density (£ s ,o = fpg ■ £0) which 
determines the final core mass. Thus it is important to keep in 
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13.5 


»I[Fe/H] 
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13.31 


13.17 
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-0.03 


C M 


9.33 x 10' 2 


1.21 x 10 13 



mind that the metallicity as it is given here serves as a proxy for 
the core mass of the object and only in our fiducial model at a 
distance of 5.2 AU it is correct to infer the M 50 value from the 
metallicity. 

The results for the fitting parameters can be found in Tab. 
4. 

Looking at the helium dependency of M50, one finds that 
the slope of both fitted lines (i.e. for hot and cold start) agrees 
quite well and approximately coincides with the values from 
SBM11, who found m Y = -20. Under variation of [D/H] and 
[Fe/H], the cold start objects react somewhat stronger than the 
hot start objects, resulting in a steeper slope. The slope for the 
hot start objects with different deuterium abundances [D/H], 
again coincides quite well with the results of SBM1 1 who found 
'«[D/H] = -1.6 X 10 4 . 

The dependence of M50 SBM11 found for variations of 
the metallicity [Fe/H] is different. This is, however, not surpris- 
ing, as the metallicity is determining only the opacity in their 
results, whereas it determines only the core mass in our results. 
It is clear that it is not expected that those two different effects 
should have the same consequences for deuterium burning. 

The most obvious difference between cold and hot start 
objects is found in the M50 sensitivity to different M max values. 
As described in Section 4.2, the deuterium burning process (in 
terms of how much deuterium is burned) is rather insensitive to 
M mal for the hot start case (an explanation for the remaining 
sensitivity will be given below). In the cold start case, however, 
the objects are much more sensitive to changes in M max . Here, 
M50 shows a variation that is about 3 times larger. In absolute 
terms, however, the variation is very small: An increase of 
M mdx by a factor of 100 (which should cover the range of M max 
expected to occur in a self-gravitationally stable accretion disk) 



only changes M50 by 



0.3? 



=2.6 



The reason for this was 



already given above: As objects with higher M max accrete a 
larger part of their mass at larger radii than the objects with 
lower M max (because due to the high accretion rate the mass 
is added at times where the objects have not yet contracted as 
much as in the low M max cases), they start their post-formation 
evolution with a higher specific entropy and have to release a 
higher amount of remaining gravitational potential energy via 
contraction (contributing to Li nt ). As outlined in Sect. 5.2 this 
makes them hotter and leads to a higher burning rate. As there is 
no entropy reducing shock in the hot start scenario, or, in other 
words, as the entire released gravitational potential energy must 
be processed out of the object by Lj nt , this effect is not observed 
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Fig. 12. Masses at which 50 % of the deuterium has been burned by 10 years. Upper row, left: varied helium fraction. Right: varied 
metallicities. Lower row, left: varied maximum gas accretion rates. Right: varied deuterium abundances. The red dashed line shows 
the fitted functions for the hot start objects (lower line in all plots), while the fits for the cold start objects are plotted in a blue solid 
line (upper line in all plots). The blue diamonds show the positions of M 5 q in the cold start runs, while the red semi-filled circles 
show the positions of M50 for the hot start runs. 



for hot start objects. 

Different maximum gas accretion rates will also lead to 
slightly different core masses in our simulations. This is a nat- 
ural consequence of the fact that the object's capture radius for 
planetesimals, and thus it's core accretion rate, depends on the 
object's structure (see Mordasini et al. 2012), and thus on M max . 
We therefore want to check whether the dependence of M50 on 
M mdx is not simply a consequence of the different associated 
core masses. However, for the cold start objects, a varying core 
mass for varying M max can be ruled out as the (main) reason for 
the sensitivity of the deuterium burning process to M max : 

In the transition region between no deuterium burned at 
all and complete deuterium fusion, all objects of one given M max 
or dust-to-gas ratio, but of different total mass, have nearly the 
same core mass. Thus we are able to plot M50 against the core 
mass for the different dust-to-gas ratios and for the different 
M max . M50 as a function of the core mass for hot and cold start 
simulations is shown in Fig. 13. If the varying amount of burned 
deuterium in the cold start case for different M max was caused 
by different core masses, then the M max simulations should lie 
on the same line as the dust-to-gas ratio simulations (the fpg 
runs). Fig. 13 shows that this is not the case, which implies 
that indeed the accretion at larger radii is important for objects 



with larger M max . For the hot start objects, in contrast, we 
find that the dependence of M50 on M max can be explained by 
different core masses alone, as the M max runs in Fig. 13 lie on 
the same line as the dust-to-gas ratio runs. It is thus not the case 
that a different M max directly influences the deuterium burning 
process in the hot start case, but it is a second order effect that 
different M max result in slightly different core masses. In order 
to get the deuterium burning mass limit M50 as a function of the 
core mass, we fit lines through the cold start runs with varying 
dust-to-gas ratio and varying M max , as well as through the hot 
start dust-to-gas ratio runs. 

M 50 /Mj 

— w dust,cold ' Moore/Me + ^dust,cold 

M 50 /Mj = mucoid • M core /M e + (29) 

"dust.hot 

This yielded the fitting parameters given in Tab. 5. 

We find that M50 is approximately equally sensitive to 
the core masses for both the cold and the hot start scenario: In 
the core mass range from ~ 30 M m to ~ 100 M m , M50 decreases 
by ~ 0.6 Mj. The y-intercept of the fitted lines furthermore 
indicates that for objects without a core M 5 o is 13.58 M s for 
a cold start and 13.42 Mj for a hot start in our fiducial model. 
However, this result should be treated cautiously, as an object 
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Fig. 13. M50 as a function of the core mass for different dust-to-gas ratios and different M max . The left panel shows the results for a 
cold start, while the the right panel is for a hot start. The black stars correspond to the runs with different M max , while the diamonds 
correspond to the runs with different dust-to-gas ratios. 



Table 5. Fitting parameters for M50 as a function of the core 
mass 



Parameter 


Value 


m dust,cold 


-9.77 x 10~ J 


^dust,cold 


13.58 


m Ai,cold 


-8.07 x 10- 2 


^Af.cold 


15.89 


m dust,hot 


-7.91 x 10~ j 


^dust,hot 


13.42 



formed via the core accretion scenario harbors a solid core at 
all times during its formation and subsequent evolution in our 
model (no core erosion) and thus we cannot further test this 
result. The second reason is that the mass limit M50 for a core 
mass of M m is an extrapolation of our data. 



6. Distinguishability of cold and hot start objects 

As both cold and hot start objects will eventually forget their 
initial conditions their evolutionary sequences will converge 
after a certain time, making an observational identification of 
the accretion mode (i.e. cold or hot start) difficult. The timescale 
after which the luminosity of a 10 Mj cold start object converges 
to the luminosity of a 10 M s hot start object is in the order of 10 8 
to 10 9 years (Marley et al. 2007). This timescale will increase 
with mass, which can be understood by looking at the specific 
entropies of the objects. For hot start objects the specific entropy 
directly after formation is increasing with mass, whereas the 
specific entropy of cold start objects is decreasing with mass. As 
the objects are nearly fully convective the specific entropy can 
be used as a proxy for the thermodynamic state of the objects. 
Thus one finds that the difference between the thermodynamic 
states of cold and hot start objects after formation increases with 
mass and thus shifts the convergence to later times. In Section 
5.1.4 we have shown that deuterium burning objects of the same 
mass would attain the same thermodynamic state if one holds 
their deuterium abundance artificially constant (the hypothetical 
"deuterium main sequence"), independent from whether they 



formed hot or cold. With non-constant deuterium abundance 
we found that hot and cold start objects will try to attain the 
"main-sequence"-state, even though they are unable to reach 
it fully. As both hot and cold start objects tend to attain the 
same "main-sequence" state it is clear that their thermodynamic 
properties will start to converge. After having burned up all 
their deuterium they will then continue the normal evolution, 
governed by contraction and cooling. This tendency of both hot 
and cold start objects to attain the same state due to deuterium 
burning should thus have an influence on the convergence 
timescale between cold and hot start objects. This is indeed 
what we find. 

First we define the timescale Tl as the timescale after which 
the ratio Lhot/^coid, where Lhot and L co id stand for the total 
luminosity of a hot and cold start object, has reached a value of 
1.1 as a function of mass. We considered two cases: For one run 
deuterium burning was disabled artificially, while for the other 
run our normal model was used, including deuterium burning. 
In our fiducial model, we found that for masses between 14 Mj 
to 21 Mj tl increases monotonously as a function of mass 
for disabled deuterium burning, going from 4 x 10 years to 
1.1 x 10 9 years. For deuterium burning objects it attains a 
constant value of approximately 2 x 10 8 years in the same 
mass range. We thus find that deuterium burning shortens the 
convergence timescale of massive cold and hot start objects. 
Nevertheless, one should keep in mind that the exact numerical 
values of 17, depend on the assumptions (e.g. the structure of 
the accretion shock) and initial conditions used for a model. 
However, the qualitative behavior should stay the same. 

7. Summary 

We have studied deuterium burning in objects with masses be- 
tween 10 to 30 Mj that form via the core accretion mechanism. 
Such high masses could be reached if gap formation does not act 
as a limit for gas accretion (Kley & Dirksen 2006), which might 
be supported by the observational finding that the planetary mass 
function seems relatively continuous up to a mass of ~ 25 Mj 
(Sahlmann et al. 2011). Our simulations self-consistently link 
the formation and the evolution phase. We have considered two 
limiting cases for the property of the accretion shock occurring 
during the gas runaway accretion phase. First, that it is com- 
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pletely radiatively inefficient, so that no liberated gravitational 
potential energy is radiated away. Instead, high entropy mate- 
rial is incorporated into the object, leading to a so-called hot 
start. We find that the properties of such objects forming grad- 
ually with the core accretion but without radiative losses at the 
shock are the same as found in the classical hot start simulations 
where one starts with an already fully formed object with some 
specified high entropy, usually thought to have formed via direct 
collapse. This shows that the structure of the shock is as deci- 
sive to determine the thermodynamic properties of an object as 
it's formation mechanism. Second, we assumed that the shock is 
fully radiatively efficient so that all potential energy liberated at 
the shock is radiated away. This leads to low entropy, so-called 
cold start objects. We have then studied the following topics: 

- Comparison with Burrows et al. (1997) and Spiegel et al. 
(2011): As seen in Section 4, we were able to cover the re- 
sults of Spiegel et al. (2011) for deuterium burning in clas- 
sical hot start objects, who started at the final mass at an ar- 
bitrary entropy without considering the formation. The radii 
at which the objects get partially stabilized against contrac- 
tion as well as the peak values of Lo/L tot and the tempera- 
tures agree very well. The only difference is that our objects 
in the transition zone between non deuterium burning ob- 
jects and deuterium burning objects seem to be able to main- 
tain a higher burning rate over a longer period of time, en- 
abling them to burn more deuterium in total. This is probably 
due to different opacities (leading to different cooling tracks) 
and the high temperature dependence of deuterium burning. 
Also, looking at the sensitivity of the minimum mass limit of 
deuterium burning in hot start objects, we find approximately 
the same slopes for M 5 o (which is the mass at which 50 % 
of the deuterium has been burned) as Spiegel et al. (2011) 
when considering variations of the helium fraction Y and the 
deuterium abundance [D/H] (see Section 5.2). In the hot start 
case we find a mass limit of M50 = 13.15 Mj for our fiducial 
model (Y = 0.25, [D/H] = 2 x 10~ 5 ). 

- Behavior of deuterium burning in hot start objects: In 
addition to the fact that we could reconfirm the results 
of Spiegel et al. (2011) we were furthermore able to in- 
vestigate the sensitivity of deuterium burning to quantities 
linked to the underlying core accretion formation process. 
We found that deuterium burning in hot start objects is vir- 
tually independent of the maximum allowed gas accretion 
rate M max and that a higher core mass (determined by the 
dust-to-gas ratio, i.e. the metallicity of the protoplanetary 
disk) yields higher central temperatures which will increase 
the deuterium burning rate. Using our fiducial model, we 
find that M50 = 13.15Mj for a core mass of 35 M e and 
M50 = 12.62Mj for a core mass of 101 M B . Extrapolation 
to M core = yields M 50 = 13.42 Mj. 

- Behavior of deuterium burning in cold start objects: As 
cold start objects accrete most of their mass at small radii 
during the runaway phase, they will be below their (hypo- 
thetical) deuterium main sequence radius when their central 
density and temperature become high enough to burn deu- 
terium. The thermostatic nature of deuterium fusion will then 
cause the objects to expand. Due to the concurrent decrease 
of deuterium they will not, however, be able to fully reach the 
hypothetical deuterium main sequence radius. This leads to 
a re-inflation of the planets to about 1.4 Rj and 2.3 Rj for for 
a 14 Mj and 23 M s object, respectively. Objects with masses 
higher than approximately 23 Mj burn almost all their deu- 
terium during the formation phase (for M max = 10~ 2 M e /yr), 



as higher masses lead to higher central temperatures and 
thus higher deuterium burning rates. Objects with masses 
< 20 Mj burn most of their deuterium after the formation 
with M max = 10~ 2 M ffl /yr is complete. 

- Sensitivity of the deuterium burning mass limit: The sen- 
sitivity of the minimum mass limit for deuterium burning for 
cold start objects, M50, is comparable to the hot start objects 
when considering variations in the helium mass fraction Y, 
the deuterium abundance [D/H] and the core mass. The de- 
pendence is, however, somewhat steeper for cold start ob- 
jects. For an increase of Y = 0.22 to Y = 0.3, M50 decreases 
from 13.9 M } to 12.3 Mj. For an increase of [D/H] from 10" 5 
to 4 x 10 ~ 5 , M 50 decreases from 13.6 M, to 12.9 Mj. When 
varying the maximum allowed gas accretion rate M max , the 
deuterium burning in cold start objects reacts stronger than 
in the hot start case. The reason for this is that the cold start 
objects with higher M max will accrete mass in the runaway 
phase at larger radii than those with smaller M max . Accretion 
onto a larger object means that the accretion shock is weaker, 
leading to less radiative losses. The object therefore has a 
higher specific entropy and, quite counterintuitively, burns 
deuterium more intensively (see Sect. 5.2 and appendix A 
and B). Similar as for the hot start, we find (using our fidu- 
cial model) that M50 = 13.3Mj for a core mass of 32 M ffi and 
M50 = 12.66Mj for a core mass of 97 M @ . Extrapolation to 
M C0Ie = yields M 50 = 13.58 M s . 

- Evolution of the luminosity: The evolution of the lumi- 
nosity (as well as the radius and effective temperature) of 
hot start objects is characterized by a monotonic decrease in 
time, once the formation process is over. Cold start objects in 
contrast begin after formation with a luminosity one to two 
orders of magnitude lower, but then get brighter, due to deu- 
terium burning. The luminosity of hot and cold start objects 
more massive than ~ 15 Mj converges after ~ 200 Myrs, due 
to deuterium burning. 

In summary, we find that the transition between deuterium burn- 
ing and non-deuterium burning objects lies at approximately 
13 Mj also for objects forming via the core accretion scenario 
(see, e.g., Fig. 12), independent of whether they might have 
formed with a hot or a cold start. The exact behavior of deu- 
terium burning, however, depends on the initial conditions such 
as the helium mass fraction, the deuterium abundance and the 
metallicity. Furthermore the exact thermodynamics of the accre- 
tion shock (cold start or hot start) and the maximum allowed 
mass accretion rate during the runaway formation phase play a 
role. Nevertheless the transition from no deuterium burning at 
all to complete fusion of all deuterium nuclei always goes from 
~ 10 Mj to ~ 15 Mj for a realistic range of parameters (see Fig. 
10). 

8. Conclusion 

Deuterium burning so far has been studied in detail only in 
objects which have not formed via the core accretion scenario 
(although Baraffe et al. 2008 found that deuterium burning 
in objects harboring a solid core is possible). It is, however, 
important to consider deuterium burning as well in objects with 
a core which have formed via the core accretion scenario, as 
it cannot be excluded that objects forming in a protoplanetary 
disk via core accretion can grow beyond the required mass of 
~ 13 Mj (Kley & Dirksen 2006). This scenario is particularly 
interesting as objects that form via the core accretion scenario 
in a disk are usually named "planets". This contradicts the 
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planet definition on deuterium burning alone, leading to an 
interesting class of deuterium burning planets (Baraffe et al. 
2008). Therefore one should consider the possibility of objects 
with masses in the Brown Dwarf domain that have a solid core, 
at least initially, as the exact fate of the core under the extreme 
conditions in objects as heavy as considered here is uncertain 
and our results are in the limit that the core is not eroded. 

Considering a core accretion formation of the objects in 
the hot start scenario, we find, when comparing to the results 
of Spiegel et al. (2011) and Burrows et al. (1997), that such 
objects with a core burn their deuterium in approximately the 
same way as the objects without a core. The sensitivity of the 
burning process to changes in the parameters like Y and [D/H] 
is approximately also the same (neglecting those parameters 
which are characteristic to the core accretion scenario). 

It is important to note, however, that if an object forms 
via core accretion and a low initial entropy (cold start) this will 
considerably alter the way in which the deuterium is burned 
and how the objects react to it. The biggest difference between 
the two formation modes is that cold start objects will expand 
to radii somewhat smaller than the hypothetical deuterium 
main sequence radius (which would be attained with a constant 
deuterium abundance). The most massive cold start object 
considered here (23 Mj) reinflates to a radius larger than 2 Rj, 
coming from a radius of just over 1 Rj. For hot start objects no 
reinflation occurs. Future direct imaging observations of young 
objects, if possible coupled with dynamical mass estimates may 
allow to observationally test these findings, and to constrain the 
formation mode of objects at the transition between planets and 
Brown Dwarfs. 
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Fig. A.l. Internal energy of the envelope obtained from the sim- 
ulation (red dotted line), the internal energy obtained by using 
the potential energy and the pressure from the simulation in Eq. 
A. 5 for the two limiting y (black solid line) and the internal en- 
ergy obtained by using the fully analytical approximation from 
Eq. A. 9 and y - 1.2857 (blue dot-dashed line). The gray area 
shows where one would expect the internal energy (red dotted 
line) to be. 

Appendix A: Envelope cooling due to the core 

The goal of this approximation is to show that even an ideal gas in the envelope 
will eventually cool in the presence of a solid core. From the usual analysis 
without a core one obtains that r cent oc l/R (see Eq. 23). In order to show this 
we make a virial analysis. Assuming virial equilibrium of the envelope we can 
start with the equation of hydrostatic equilibrium 



In addition to the results of Spiegel et al. (2011), who 
found that deuterium burning is depending on the metallicity 
of the objects via the opacity, it should be kept in mind that 
the metallicity also determines the core mass (if the local gas 
surface density in the protoplanetary disk is known). A higher 
core mass leads to a higher burning rate. Thus this should add 
up to the opacity related deuterium burning enhancement effect 
of the metallicity. 

Finally we saw that the consideration of a core accretion 
scenario also yields that the border between deuterium burning 
and non deuterium burning objects remains at approximately 13 
Mj . To model the exact deuterium burning behavior, however, it 
is important to know certain characteristics of a system, such as 
the helium abundance, the solid surface density, the deuterium 
abundance and an upper limit for the gas runaway accretion 
rate. 

Numerical data showing the formation and subsequent 
evolution of cold and hot start objects in our fiducial model can 
be found at http://www.mpia.de/homes/mordasini/Site7.html. 



Acknowledgements. We thank Hubert Klahr, Yann Alibert and Kai-Martin 
Dittkrist for helpful discussions. We especially thank David Spiegel and Adam 
Burrows for dedicated comparison calculations concerning the expansion phase 
of cold start objects. Finally, we thank our referee Gilles Chabrier for questions 
which helped us to gain a deeper understanding of the governing physical mech- 
anisms. Christoph Mordasini acknowledges the support as an Alexander von 
Humboldt fellow. 



= — VP - V(j> 
P 



(A.l) 



where P is the pressure, p the density and (/> the gravitational potential. We as- 
sume the gas in the envelope to be ideal. Furthermore we assume that the core 
density is constant. In this idealized analysis, one writes 
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We integrate over the envelope only as we are interested in the envelope only. 
The core is merely treated as an external factor; this also means that the second 
integral on the RHS of the above equation contains the potential energy of the 
envelope in the gravitational field of the core, but not the potential energy of the 
core itself. Using dm/p = dV the first integral on the RHS of Eq. A.2 transforms 
to 

r-VP— = -<£ PrdS +3(y- 1) f pudV , (A.3) 

P J env J env 

-3P s V pl +3i>cV c fin 

where dS is an infinitesimal vectorial surface element of the envelope, y is the 
adiabatic index and u is the specific energy density of the gas. P s is the surface 
pressure which is negligible in our problem, V v \ is the object's volume, P c is 
the gas pressure on top of the core, V c is the core volume and U mt denotes the 
internal energy of the envelope. Furthermore one can show that the second term 
on the RHS of Eq. A.2 is equal to the gravitational potential energy W of the 
envelope with 

W = Wenv-core + W cn v-env (A.4) 

where Wenv-core ' s tne potential energy of the envelope in the gravitational field 
of the core and W C nv-env is the potential energy of the envelope due to its own 
gravitational field. Thus we finally have 
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We furthermore found that W cnv _ clw is reasonably well described by 
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Wenv-core can be estimated by 

" env— ei 
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To approximate the core pressure we use the equation of hydrostatic equilibrium 
again and neglect all terms linear in the core mass. One then finds 
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Again using that P s is negligible one finally gets by using Eq. A. 5, A. 6 and A. 
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With this approach one finds that the envelope starts to cool if the R p \ becomes 
smaller than 2.83 R c if M c — > and somewhat later for non negligible core 
masses. The physical reason for this cooling process is that the excess energy is 
transferred to the core which has to withstand an evergrowing external pressure. 
In order to test this approximation we carried out a simulation forming a 1 Mj 
planet, using a constant core density as well as an ideal H9 EOS for the envelope. 
y then varies between 1 .4 near the surface and 1 .2857 at the center. The core mass 
was 35 M e . In Fig. A. 1 one can see the internal energy of the envelope obtained 
from the simulation, the internal energy obtained by using the potential energy 
and the pressure from the simulation in Eq. A. 5 and the internal energy obtained 
by using the fully analytical approximation from Eq. A.9 and y = 1.2857. One 
sees that the analytical approximations recover the shape of the internal energy 
curve obtained from the simulation quite well. Furthermore the simulation result 
is bracketed by the two results obtained from Eq. A. 5 using y = 1.4 and y = 
1.2857. 



Appendix B: Envelope cooling due to the core and 
the electron degeneracy 

It is well known that the ions in fully degenerate objects, such as White Dwarfs, 
will cool as they contract. The excess energy is put into the degenerate electrons. 
In the objects considered in this paper we have partially degenerate envelopes. 
However, in order to make an analytical analysis simpler, we will consider a fully 
degenerate envelope. In combination with the effect that the presence of the core 
can already cool an ideal gas envelope we expect a cooling behaviour for the 
degenerate envelope as well. We use the simplification that the internal energy 
of the electrons is approximately 



1 N.E F 



(B.l) 



where Ef is the Fermi energy and N e the number of electrons. Next we use that 
Ef oc p 2 F and p F oc n, , where pf is the Fermi momentum and n e the electron 

number density. As n c oc (r 3 — R^\ one finally obtains that 
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Here we made the approximation that R?^ a j{3 - which is justified in the 
case considered here as degenerate envelopes do not contract as much as the 
ideal ones do. Furthermore it holds for the potential energy of the envelope 
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which yields, together with Eq. B.2: 
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We use that the ions can be considered as an ideal gas. Then it holds for both the 
ions and the degenerate electrons that P = 2/3p«. One finds again from a virial 
analysis that one can write then 
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Here we neglected W cnv _ colc this time since for the objects we are ultimately 
interested in, it holds that M em s> M c . As we consider a fully degenerate enve- 
lope the internal energy of the electrons should be much larger than the internal 
energy of the ions: 

l/ intc » U- mtl . (B.6) 



Thus one finds from Eq. B.5: 
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From Eq. A. 6 and Eq. A.8 we know that 

'iPcVc = -W em - em K R (R p \,R c ) , 

with 

R P i - «c «c 

For the change of the internal energy dU mt = dll mtt . + dU mtl we now find 

dUinu + dt/tot, = -^env-env " ^d(3P c V c ) . (B.10) 
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Using Eq. B.4 and Eq. B.8 then yields 
dU m 1 ' 



, K R \dW em ^ m + -d(W sm - aw K R ). (B.ll) 



As the only free parameter for a given object is R p \ (i.e. M em and R c are fixed) 
one finds the above equation to be equal to 
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For the cooling of the ions one must have 
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which means that the internal energy of the ions decreases as the objects shrink. 
In the left panel of Fig. B.l one can see that this is indeed the case, i.e. we expect 
the degenerate envelope around the solid core to cool. For a partly degenerate 
envelope however, if the objects radius is much larger than the cooling radius 
derived for ideal gas (2.83 R c ), the object should heat up during contraction. The 
reason for this is that the degeneracy in the envelope is not yet strong enough to 
result in ion cooling. In the right panel of Fig. B.l one sees the central temper- 
ature (i.e. the temperature in the envelope above the core) of a 10 Mj hot start 
planet as a function of the radius for the whole post-formation phase (starting at 
time to). The envelope EOS is the EOS of Saumon et al. (1995), and the core 
density was not kept constant artificially. One can see that the envelope heats up 
until it has reached a radius of approximately 1 .7 Rj and then cools down. Non 
deuterium burning cold start objects will always cool in the post formation phase 
as they are already strongly degenerate. 
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